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Abstract 

The properties of large bipolarons in two and three dimensions are investigated by averaging over the 
relative wavefunction of the two electrons and using the Lee-Low-Pines-Huybrechts variational method. 
We obtain the ground-state (GS) and excited-state energies of the Frohlich bipolaron for the whole range 
of electron-phonon coupling constants. Furthermore, we calculate the energies of the first relaxed excited 
state (RES) and Franck-Condon (FC) excited state of the bipolaron. Compared with the FC state, the first 
RES has a lower energy. Our results for the GS and RES energies are lower than those obtained before by 
the Landau-Pekar method in the whole coupling regime. 
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I. INTRODUCTION 



Polaronic anc bipolaronic effects are found in many materials such as transition metal oxides 
111, polymers J2], and superconducting materials fl. Since these materials possess low dimen- 
sional structures, the effective electron-phonon coupling is very strong owing to the confinement 
of their low dimensional nature. With this feature, many intriguing phenomena related to polarons 
and bipolarons can occur in these materials. Electrons in a polar or ionic crystal are influenced 
by the phonon field giving rise to polaronic effects When the attractive interaction of 

the electron-phonon coupling is strong enough to overcome the Coulomb repulsion between two 
electrons, the two electrons plus their associated phonon field can form a bound state. Such a com- 
posite quasi-particle is referred to as a bipolaron. If the bipolaron is extended over many lattice 
cells, it is called large bipolaron. In contrast, if the bipolaron is localized at one lattice site, it is 
called small bipolaron [6]. 

Bipolarons have been extensively discussed both for fundamental theoretical reasons and for 
their importance in semiconductor materials. With recent advances in the creation of nanocrystals 
and semiconductor nanostructures, strong electron-phonon coupling is realized due to quantum 



confinement effects 
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1311 . and the properties of bipolarons in low dimen- 



sional systems are of growing interest. This problem is also relevant to the proposal of the bipo 



laronic mechanism for electron pairing in the Cu0 2 plane in high T c cuprates |3J, |6J. Recently, 
Schoenes and co-workers ful reported that some superconductors such as Y Ba 2 Cu^0 8 have 
strong electron-phonon coupling, which also boosts the theoretical research on bipolarons. 
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Though the bipolaron ground state has been extensively discussed in the past decade H15LI16L 
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2211 . the excited states have, so far, been studied much less. At the same time, 



knowledge of the excited states of the bipolaron is related to electron transport, photolumines- 

JOES . Several works have been devoted to excited states of Frohlich 
, 22D. Huybrechts I23I developed a Lee-Low-Pines variational method 



cence and photoemission 



bipolarons ||23 
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and studied the ground and excited states of a single optical polaron. Smondyrev et al. [24] cal- 
culated the energv spectra of the one-dimensional bipolaron in the strong-coupling limit. More 
recently, Sahoo [25] developed the Landau-Pekar variational method to get the ground and first 
excited states of the Frohlich bipolaron in a multidimensional ionic crystal in the strong-coupling 
limit. 

In this paper we adopt the Frohlich Hamiltonian for the large bipolaron which consists of 
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the kinetic energies of two electrons, their interactions with the phonon field, and the screened 
Coulomb repulsion between them. We extend the Huybrechts variational approach (LLP-H) to the 
analysis of the bipolaron ground state in Sec. II. For comparison, different wavefunctions are used 
to calculate the ground- state energy of the bipolaron in Sec. Ill, and the best wavefunction of the 
relative motion is obtained. In Sec. IV, after averaging over the relative motion, the energies of 
the first relaxed excited state and Franck-Condon excited state of the bipolaron in two and three 
dimensions are obtained for the whole coupling parameter range, and discussions of these two 
types of excited states are given. Finally, comparison of our results with previously published data 
is given in Sec. V. Our conclusions are presented in Sec. VI. 

II. CALCULATION OF THE GROUND-STATE ENERGY 

The Frohlich Hamiltonian for the polaron in N dimensions (ND) has been derived by Peeters 
et al J^sl- Accordingly, the Hamiltonian describing a system of two electrons interacting with a 
longitudinal optical (LO) phonon field may be written as 



where all vectors are iV-dimensional (N = 2, 3), r-j (p^) is the position (momentum) operator of the 
ith electron (i = 1,2) and m is the effective electron band mass in the parabolic approximation. 
a£ and a k are respectively the creation and annihilation operators of the LO phonons with the 
wave vector k. Here we should mention that the impurity-phonon interactions have already been 
eliminated so that we assume c<j k = uj LO - For convenience in the following, units are taken where 
m — h — u L o — 1, and the Hamiltonian for the ND Frohlich bipolaron has the form: 
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In this way, energies are in units of huj LO and lengths in units of ^Jh/imuio) in this article. The 



interaction coefficient is 
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where Vn is the ND crystal volume, a is the dimensionless electron-phonon coupling constant 

1 e 2 / J_ _ J_\ /2m^o\ 1/2 ^ (4) 
/iculo 2 V^oo So) V ^ / 
and (e ) is the high-frequency (static) dielectric constant of the medium. U(r) = U/r is the 
Coulomb interaction potential between the two electrons where the nonscreened electron Coulomb 
repulsion strength is given by U — e 2 /e OQ , which may be rewritten as 

TT ^ a f<\ 
U = - , T] = . (5) 

1-r) e 

Since the bipolaron is a composite particle, it is convenient to introduce center-of-mass and 
relative coordinates and momenta, R = (ri + r 2 )/2, P = pi + p 2 , r = ri — r 2 , p = (pi — p 2 )/2, 
in which the Hamiltonian can be rewritten as 

H= ^ + 2^cos^V^ k a k e ikR + / i .c.) 
k V - / 



u 



+ J]a+a k + p 2 + -. (6) 

k 

Up to now, an exact analytical solution for the Hamiltonian (6) is not known. However, we 
can make progress by averaging (6) over the relative wavefunction 0(r), which yields an effective 
Hamiltonian for the center of mass motion, 

H °ff= ^ + E(£ka k e tk - R + /i.c.) 

k 

+ ^a+a k + £ r . (7) 

k 

Here B k = 2V k (cos( k f )), E r = (p 2 + U/r), with (• ■ ■) denoting an averaging over the wavefunc- 
tion (f>(r). The new Hamiltonian (7) would be equivalent to (6) if the true wavefunction would be 
known which of course is not the case. Therefore, we will use a variational approach, two specific 
choices will be discussed in section III. 

Note that the effective Hamiltonian (7) corresponding to the center-of-mass motion is in essence 
equivalent to a single-polaron Hamiltonian. The differences between the Hamiltonian (7) and the 
usual Hamiltonian of a single polaron are the following: (i) the energy is shifted by the average 
of the energy of the relative motion E r , and (ii) the electron-phonon interaction coefficient £? k is 
renormalized. 

In this paper, we will follow the modified Lee-Low-Pines variational method proposed by Huy- 
brechts 12311 for the polaron problem. This will allow us to obtain very accurate bipolaron energies 



for all coupling regimes and calculate the ground and excited- state energies of the Frohlich bipo- 
laron. Performing the unitary transformation JJ\ 



Ui = exp ^— ia k • Ra^kj , 



the Hamiltonian (7) can be transformed into 



H. 



f ^(P-a^ka+a k 



+ $^{^k a k exp[-i(l - a)k • R] + h.c.} 

k 



(8) 



(9) 



where a is a parameter. In the limit a — > 0, the present calculation is identical to the strong- 
coupling regime, whereas the case a — > 1 corresponds to weak electron-phonon interaction. 

Following Huybrechts, we introduce the creation and annihilation operators 5+ and bj for elec- 
trons 



Pi = W (bj + ty), 

/ 1 \ V2 



(10) 



where the index j refers to the N space directions and A is a variational parameter. Rewriting the 
Hamiltonian (9) using (10) one gets 



H ef f = H ± + H 2 + E r , 
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(12) 
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H 2 = — k ■ k a^a^a^ay,' 
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(13) 
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Performing now the second Lee-Low-Pines transformation: 



U 2 = exp 



^(/k<£ - /k«k) 



(14) 



one obtains 



a 2 k 2 



4 + /k)(«k + / k ) 



+ E ^k(4 + /k)exp 



l-a) 2 A; 2 
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xexp 



xcxp 



1 - a) 



+H 1 + -E r , 



(15) 



where the part if x of the Hamiltonian (15) contains terms of no importance for the further cal- 
culation and will be omitted. The wavefunction describes the polaron system \I/(r) including an 
electron wavefunction 0(r) and aphonon field wavefunction |/), \I/(r) = <p(r) \ f). For the ground 
state, one has 



with 



Mr) = Mr) I/) , 



bj<fro( r ) = 0, <j> (r) = cexp(-Ay^ R)), 



(16) 



(17) 



where bj and A are defined in equation (10), <po(r) corresponds to the ground-state wavefunction 
for the electron and j refers to the jth direction of N space. For the field function in the ground 
state one has 



\f) = U 2 |0> , a k |0) = 0, 



(18) 



where {7 2 is defined in equation (14) and fk in {7 2 will be obtained by minimizing the ground-state 
energy of the bipolaron, 
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BP 



0. 



(19) 



Then we can directly present the bipolaron ground-state energy E^p D as 

a 2 k 2 ^ 
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+E r , 
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1 - a) 2 k 2 
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+ h.c. 



(20) 



with 



/a- 



-5£exp[-(l-a) 2 A; 2 /(8A)] 



(21) 



1 + a 2 fc 2 /4 

Using (21) in (20) we obtain the ground- state energy of the bipolaron in iV dimensions with an 
arbitrary electron-phonon coupling constant: 



E 



BP 



4 

+E r . 



+ a 2 k 2 /A 



(22) 



The effective electron mass of ND materials which is modified by bipolaronic effect is given by 
(for details see Ref. Q) 



, i?gexp[-(l-a) 2 fc 2 /(4A)] 
(l + a 2 fc 2 /4) 3 



(23) 



To obtain explicit results for the ground-state energy and the effective mass we now compute the 
relative wavefunction of the electron pair. 



III. CHOICE OF THE BEST WAVEFUNCTION FOR RELATIVE MOTION 

In order to calculate the energy related to the relative motion E r and the coefficient B^, we 
have to choose a trial wavefunction <p(r). In the following two different functional forms will be 
chosen: (a) Coulombic type 

</4 c )(r) = CV n exp(-fir/4), (24) 
7 
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FIG. 1: 3D bipolaron ground-state energy versus a using (a) Coulombic wavefunctions (24) with n = 
0, 1, 2, 3 for r\ = 0; (b) oscillator wavefunctions (25) with n = 0, 1, 2 for 77 = (units of energy: Hujlo)- 

and (b) oscillator type 

0(°)(r)=CV 1 exp(-5r 2 /4), (25) 

where 6 and f2 are variational parameters and n — 0, 1, 2, .... Using the wavefunctions (24) or (25) 
to calculate £V and i?k, we obtain the ground-state energy of the bipolaron from equation (22). 
Comparing the results obtained with wave functions (a) and (b) allows us to determine which one 
leads to the best approximation of the relative motion. 

In Fig. 1(a), the ground-state energy E PP D of a 3D system is plotted versus a for different 
Coulombic wave functions (24) for the case n = 0, 1, 2, 3. Interestingly, the energy obtained with 
the function 0^ (n = 2) is smaller than those obtained by the others with n = 0, 1, 3. From 
the viewpoint of the variational principle, the relative motion is best described by the function 
r 2 exp(— fir/4) among the different Coulombic wave functions. 

In Fig. 1(b), we plot E°/ P D (3D) as a function of a for different oscillator wavefunctions (25) 
for n = 0, 1,2. Here, the energy obtained with (pf is lower than those by others, indicating 
that r exp(— br 2 /A) is superior to other oscillator wavefunctions in describing the relative motion. 
Furthermore, we find that Eg^P calculated with the oscillator wavefunction ^ is smaller than 
that obtained by using the best Coulombic wavefunction 0^ • So we conclude that the oscillator 
wavefunction (p^ reflects the relative motion best in 3D materials. 

Next consider the 2D case. Fig. 2(a) and Fig. 2(b) present Eg 2 p D as a function of a for Coulom- 
bic and oscillator wavefunctions, respectively. In Fig. 2(a), we find that the energy obtained by 
0^ is lower than those by (p^ and 0^- Fig- 2(b) shows that the energy obtained from 0^ is 
lower than those following from and 0^. Compared with E PP D of 0^, E PP D of (p^ is a 
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FIG. 2: 2D bipolaron ground-state energy as a function of a (a) for Coulombic wavef unctions (24) with 
n = 0, 1, 2 at j] = 0; (b) for oscillator wavefunctions (25) with n = 0, 1, 2 at 77 = (units of energy: huJLo)- 

little smaller, which demonstrates that the oscillator wavefunction (jy^ reflects the relative motion 
best in 2D systems as well. 

Finally we calculate the energy of the relative motion E r by using the best wavefunction (f)^ 
in equation (22) and obtain the ground-state energies of the bipolaron in 3D and 2D materials as 
follows 

k 2 k 4 \ 



E ° p =T {1 - a) -—J dk { 1_ 66 + 1446* J 
exp(-gk 2 ) 7b 4 j~b~ 
X l + a 2 P/4 + 12 + 3 V 2^ 

-(l-a) 2 -2V2a [°° dk ( 1 - ^ + 
2 V ' Jo \ 4b 64b 2 

exp(-gk 2 ) b V2to 
X l + a 2 fc 2 /4 + 2 + 4 ' 



(26) 



ciO,2D 
^BP 



(27) 



with 



^=46 



4A 



(28) 



and a, 6 and A being usual variational parameters. Comparing our bipolaron ground-state energy 
(Egp D ) with the ground-state energy (Ep' ND ) of a single polaron 112 8h . we find that the bipolaron 
is stable for a > a c = 8.1(3.6) in 3D (2D). It is interesting to note that our results are very close 
to the estimate a c — 8.1 in 3D and a c = 3.5 in 2D in Ref. Ill 811 by Luczak et ah, who obtained the 
bipolaron ground-state energy by a totally different variational scheme. 

Furthermore, the effective mass and root mean square separation r 12 of two electrons can be 
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FIG. 3: Effective electron mass in two and three dimensions as a function of the coupling constant at 
T] = 0.08 (units of bipolaron effective mass: free electron mass m e ). 



calculated as 



m 



*3D 



2 + 



x 



A^2a 
3tt 

k 2 exp(—gk 2 ) 
(l + a 2 A; 2 /4) 3 ' 



dk 1 



6b + 1446 2 " 



(29) 



m 



*2D 



2 + V2aj^ dk(^l 

k 2 exp(—gk 2 ) 
X (l + a 2 £; 2 /4) 3 ' 



P k^_ 
46 + W 



(30) 



,.3D 
12 



1/2 



(31) 



„2D 
12 



1/2 



(32) 



The values of the variational parameters a, b and A which are used in (26-32) are presented in 
tables 1. The effective mass as a function of the coupling constant is shown in Fig. 3. Due to 
the interaction of electron and phonon field, the effective mass greatly increases with the coupling 
constant. It also can be seen from Fig. 3 that the effective mass in 2D materials is larger than that 
in 3D materials, which indicates that the electron-phonon interaction is increased as the dimension 
of the material is reduced. Similarly, in one-dimensional materials H2411 . the bipolaron has a much 
larger effective mass than our results of 2D materials. 
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TABLE I: The values of the variational parameters A, a and b which are used in the GS and FC state energies 
for 3D and 2D materials. 



a A 3D a 3D b 3D X 2D a 2D b 2D 



I. 0.628 0.802 0.7297 1.58 0.635 1.795 
2.0 1.52 0.638 1.659 4.26 0.376 4.743 
3.0 2.60 0.507 2.715 8.51 0.217 9.396 
4.0 3.93 0.398 3.990 14.6 0.132 16.11 
5.0 5.58 0.310 5.574 22.6 0.0871 24.88 
6.0 7.61 0.241 7.530 32.4 0.0612 35.68 
7.0 10.1 0.188 9.889 44.1 0.0452 48.47 
8.0 12.9 0.150 12.66 57.5 0.0347 63.25 
9.0 16.2 0.121 15.84 72.7 0.0274 80.02 
10.0 19.9 0.0991 19.42 89.8 0.0222 98.76 

II. 23.96 0.0826 23.40 108.6 0.0184 119.5 
12.0 28.4 0.0698 27.77 129.2 0.0155 142.2 
13.0 33.3 0.0597 32.53 151.7 0.0132 166.8 
14.0 38.6 0.0516 37.68 175.9 0.0114 193.5 
15.0 44.3 0.0450 43.22 201.9 0.0099 222.1 



IV. INTERNAL EXCITED STATES 

Following the definition given by Devreese fl, we compute the relaxed excited state (RES) and 
Franck-Condon excited state (FC) energies of the bipolaron. The RES is created if the electron 
in the polaron is excited while the lattice readapts to the new electronic configuration. One can 
imagine that the electron goes from a Is- to a 2p-state, while the lattice polarization in the final 
state is adapted to the 2p-state of the electron. If, on the contrary, the lattice corresponds to the 
electron ground state, while the electron is excited, one speaks of a FC state. 



Here we use a method of Huybrechts 112311 who computed the FC state energy 13011 of a single 
polaron. In this section, we will develop Huybrechts method to calculate the RES energy of the 
bipolaron by adjusting a, b, A and fa to the relaxed excited state. Finally the FC state and RES 
energies we obtain are fit for the whole electron-phonon coupling regime. 
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The wavefunction for the excited polaron system is as follows 



with 



*i(r) = 0i(r) |/) , 



b+Mr) = Mr), Mr) = cR 3 exp(-A^ R 2 ). 



(33) 



(34) 



In the FC model, the electrons can be excited to a higher level in the same potential well built 
up by the field of virtual phonons in the ground state so that A has the same value as that in Tab. 
1. For the RES, A is different from that in Tab. 1 and it can be obtained by minimizing the first 
excited-state energy. The potential well should adjust to the first excited state, Consequently, the 
lattice polarization in the final state is adapted. 

For the field wavefunction in the excited state one has 



I/) = U 2 |0> , a k |0> = 0. 



(35) 



For the FC state f k in U 2 is the same as that in equation (21). For the RES f k can be obtained by 
minimizing the excited-state energy according to 



dE 



1,ND 
BP 



= 0. 



(36) 



In order to obtain the energy of the first internal excited state of the bipolaron we calculate the 
average 



(0 | Mr) I H" | Mr) I 0). 



(37) 



After some elementary calculus, we get the excited-state energy E^p D for the bipolaron in N 
dimensions 

a 2 k 2s 

k 



hND A iVA(l-a) 2 
h BP = - H ; 



+E *£/. 



k 

x exp 



1 P 

l-(l-a) 2 — - 
V 4A iV 



1 - a) 2 k 2 
8A 



+ h.c. } + E r , 



with 



FC _ -^ex P [-(l-a) 2 P/(8A)] 



1 + a 2 A; 2 /4 



(38) 



(39) 
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for the FC state and for the RES one has 



fk 



RES _ -5£exp[-(l-a) 2 fc 2 /(8A)] 



x 



1 + a 2 k 2 /A 
1 k 2 



(40) 



(o) 



Then we calculate the energy of the relative motion by using the best relative wavefunction <p\ 
in equation (25) and inserting f£° and f^ ES into equation (38). Finally we obtain the FC state 
and RES energies of the bipolaron in 3D and 2D materials as follows 



E 



FC,3D 
BP 



A 3A . 
= 2+ T (l-«) 



4-\/2. 



x 



7T 7 V 66 14462 / l + a 2 ^ 2 / 4 
(l-o) 2 A;^ 



12A 



76 4 / 6 TT 
+ 12 + 3V^ 



(41) 
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A 
2 
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-2V2a [ 
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X 



(l-a) 2 A; 2 
8A" 



46 
6 

+ 2 + 



646 2 / l + a 2 A; 2 /4 



(42) 



BP 



A 3A 9 
= 2 + T (1 - a) 



a- 



x 



7T 



1 - 



ii. | | A: 2 A: 4 A exp(-^fc 2 ) 
' 11 66 1446 2 ) l + a 2 k 2 /4 



(1 - a) 2 k 2 
12A 



76 4 . , , 
+ 12 + 3V^ 




(43) 
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46 6 46 2 y 1 + a 2 P/4 

6 v^,, 
f 2 + — ^ 



where E^p D is expressed in units of Jtwlo- For the FC state the values of the parameters a, 6 and 



(44) 



A are taken from the calculation of the ground state energy E^p D (Tab. 1). For the RES the values 
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TABLE II: The values of the variational parameters A, a and b which are used in the RES energy for 3D and 
2D materials. 



a 


A 3D 


a 3D 


b 3D 


X 2D 


a 2D 


b 2D 


1.0 


1.13 x 10" 


- g 0.999 


0.631 1.01 x 10~ 9 


0.999 


1.227 


2.0 


3.71 x 10" 


- 9 0.999 


1 1 70 1 ; 

1.11 U ± . < 


31 x 10 -9 


0.999 


2.097 


3.0 


2.59 x 10" 


9 0.999 


1 556 


1 34 


0.717 


3.193 


4.0 


0.726 


0.792 


2 D85 


3 25 


0.515 


4.671 


5.0 


1.54 


0.651 


2 685 


5 53 


0.387 


6.478 


6.0 


2.42 


0.551 


3 337 


8 31 


0.296 


8.694 


7.0 


3.39 


0.472 


4 D65 


1 1 6 


0.230 


11.36 


8.0 


4.49 


0.406 


4 886 


1 5 5 


0.182 


14.48 


9.0 


5.74 


0.350 


5.812 


19.9 


0.147 


18.06 


10.0 


7.13 


0.302 


6.855 


24.9 


0.121 


22.09 


11.0 


8.67 


0.262 


8.020 


30.3 


0.100 


26.57 


12.0 


10.4 


0.228 


9.312 


36.4 


0.0847 31.50 


13.0 


12.2 


0.199 


10.73 


42.9 


0.0724 36.86 


14.0 


14.3 


0.175 


12.28 


50.0 


0.0625 42.67 


15.0 


16.5 


0.155 


13.95 


57.7 


0.0545 48.91 



of the parameters A, a and b are obtained by minimizing the excited-state energy E BP ' of Eq. 
(36) through numerical calculation, which are presented in tables 2. From Tabs. 1 and 2, we can 
see that A in the RES is smaller than in the ground state which means that there is weaker potenial 
well built up by the phonon field in the RES. 

In Fig. 4(a) and Fig. 4(b), the ground-state energy (E^p D ), first RES energy (E^p S ' ND ) 
and FC state energy (E BP ' ) of the bipolaron are displayed as functions of the electron-phonon 
coupling constants for 3D and 2D materials. We find that in the whole range of electron-phonon 
coupling constants, the RES energies are negative. On the other hand, the FC state energies in 3D 
and 2D materials are negative when a > 0.4 and a > 0.2, respectively. The difference in energy, 
AE ND = E B ^ D - E° B ^ D , yields the excitation energy, which is related to optical absorption of 
bipolarons in semiconductor materials J3II31I1. 

In all cases, we find that the RES energy is smaller than the FC state energy, see Fig. 4, which 
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FIG. 4: (a) Ground-state energy (E°/ p D ), first RES energy (E pp S ' 3D ) and FC state energy (E pp ' 3D ) of 
the bipolaron are displayed as functions of the electron-phonon coupling constants for 3D materials. The 

7~) 2Z) RES 2.D 

bipolaron is stable for a > a c = 8.1. (b) Ground-state energy (E B ' p ), first RES energy (E Bp ' ) 

FC ID 

and FC state energy (E Bp ' ) of the bipolaron are displayed as functions of the electron-phonon coupling 
constants for 2D materials. The bipolaron is stable for a > a^ D = 3.6. 



is well known from absorption spectrum calculations ||4J, |32|, |33|] : The absorption peak due to a 
transition from the ground state to the first relaxed excited state corresponds to the zero-phonon 
peak. In contrast, an absorption transition from the GS to the FC state is accompanied with phonon 
emission. If the bipolaron system is excited to the FC state, the lattice will relax towards the RES 
by emission of phonons. 



V. COMPARISON WITH STRONG-COUPLING CALCULATIONS 

In order to assess the heuristic value of our approach, we compare the present GS and RES 
energies with those obtained by Sahoo J25I1 . who adopted the wavefunction o °^( r ) an( ^ developed 
a Landau-Pekar variational method to get the ground- state and RES energies of the Frohlich bipo- 
laron in the strong-coupling limit. For the case of N dimensons we obtain from Eqs. (22, 25, 38, 
40) the following results when O °^ (r) is used 



4 v ) 4 v / 



E^r(0,a) = — (1 - aY + — + ^b/2U y ^ 
2v / 2«r(^i) [°°exp(-gk 2 



T(f) J l + a 2 F/4 



dk, (45) 
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4 v 4 v , r( „. 



v^F r(f) 7 l + a 2 F/4 



x[1 - L -5^- 1<tt+ 2" <46) 

with a, 6 and A being usual variational parameters, and g being defined in equation (28). If we set 
a = and use the units 2m — h — u — 1 equations (45, 46) reduce to the results of Sahoo [25]. 
To perform a comparison with our results, we define the relative deviation £ = (E^£ D (0, 0) — 
E^p D )/ |£°^ D (0, 0) |, with E°£ D (Q, 0) referring to the GS energies obtained by Sahooj25J] and 
E°^ D given by Eqs. (26,27). Also, we define 5 = ^(0, 0V-^ S '^)/ 1^^(0,0)1, 
5,7VD (0, 0) referring to the RES energies obtained by Sahoo and E*p S ' ND denoting our 
results (43, 44). 

In Fig. 5, we plot £ and S as functions of a for 3D and 2D materials at 7] = 0. £ and <5 
are positive in the whole coupling regime, indicating that our GS and RES energies are smaller 
than those of Sahoo. It is well known that the Landau-Pekar method works well only in the 
strong-coupling regime, indeed £ and 5 decrease monotonously with the increase of a. At large a, 
E®' S ^ D (0, 0) and E^p S,ND (0, 0) are very close to our results, which demonstrates the reliability of 
our approach. These numerical results demonstrate that our results which are obtained by the LLP- 

n 

H method give improved GS and RES energies compared to previous strong-coupling models[25]. 
Moreover, our results apply to the whole electron-phonon coupling regime which is due to the 
additional parameter a in equation (8). 



VI. CONCLUSIONS 



We have extended the Huybrechts variational approach (LLP-H) to the analysis of bipolarons. 
By performing an average over the wavefunction of the relative motion of the two electrons, the 
ground and first excited-state energies of the bipolaron in two and three dimensions are obtained. 
Numerical results show that the RES energy is smaller than the FC state energy. 

The energies we obtain are applicable to arbitrary values of the electron-phonon coupling con- 
stants. Our ground-state and RES energies are lower than the previously reported results from 
Landau-Pekar method iQ], which is due to the use of a more appropriate relative wave function 
and the additional parameter a in LLP-H method. Based on analytical and numerical results, we 
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h BP )l 



El? c D '(0,0) 



FIG. 5: Relative deviation of our ground-state energy and Sahoo's result (f = (eI£ D '(0,0) - 
) as a function of a for 3D and 2D materials at rj = 0; Relative deviation of our 

^RES,ND^ I t-,RE S ,N D / n A \ , f . 

i BP ) / -£/ osc (0, 0) ) as a function of a 
for 3D and 2D materials at rj = 0. 



RES energy and Sahoo's result (5 = (E^ S ' ND (0, 0) 



conclude that the best relative wavefunction is the oscillator-type function ^ both for the ground 
state and the first excited state of the bipolaron. Our results could be of relevance for high-T c 
superconductors where bipolarons are expected to play an important role [6]. 
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